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Abstract 

We  study  the  approximation  of  the  smallest  eigenvalue  of  a  Sturm-Liouville  prob¬ 
lem  in  the  classical  and  quantum  settings.  We  consider  a  univariate  Sturm-Liouville 
eigenvalue  problem  with  a  nonnegative  function  q  from  the  class  C'2([0, 1])  and  study 
the  minimal  number  n(e)  of  function  evaluations  or  queries  that  are  necessary  to  com¬ 
pute  an  gr-approximation  of  the  smallest  eigenvalue.  We  prove  that  n(e)  =  ©(M1/2)  in 
the  (deterministic)  worst  case  setting,  and  n(e)  =  0(£~2/5)  in  the  randomized  setting. 
The  quantum  setting  offers  a  polynomial  speedup  with  bit  queries  and  an  exponential 
speedup  with  power  queries.  Bit  queries  are  similar  to  the  oracle  calls  used  in  Grover’s 
algorithm  appropriately  extended  to  real  valued  functions.  Power  queries  are  used  for 
a  number  of  problems  including  phase  estimation.  They  are  obtained  by  considering 
the  propagator  of  the  discretized  system  at  a  number  of  different  time  moments.  They 
allow  us  to  use  powers  of  the  unitary  matrix  exp(^iAf),  where  M  is  an  n  x  n  matrix 
obtained  from  the  standard  discretization  of  the  Sturm-Liouville  differential  operator. 
The  quantum  implementation  of  power  queries  by  a  number  of  elementary  quantum 
gates  that  is  polylog  in  n  is  an  open  issue. 

In  particular,  we  show  how  to  compute  an  ^-approximation  with  probability  | 
using  n(e )  =  @(e_1//3)  bit  queries.  For  power  queries,  we  use  the  phase  estimation 
algorithm  as  a  basic  tool  and  present  the  algorithm  that  solves  the  problem  using  n(e)  = 
0(log£_1)  power  queries,  log2  e_1  quantum  operations,  and  |  loge^1  quantum  bits.  We 
also  prove  that  the  minimal  number  of  qubits  needed  for  this  problem  (regardless  of 
the  kind  of  queries  used)  is  at  least  roughly  ^  logs-1.  The  lower  bound  on  the  number 
of  quantum  queries  is  proven  in  [5]. 

We  derive  a  formula  that  relates  the  Sturm-Liouville  eigenvalue  problem  to  a 
weighted  integration  problem.  Many  computational  problems  may  be  recast  as  this 
weighted  integration  problem,  which  allows  us  to  solve  them  with  a  polylog  number  of 
power  queries.  Examples  include  Grover’s  search,  the  approximation  of  the  Boolean 
mean,  NP-complete  problems,  and  many  multivariate  integration  problems.  In  this 
paper  we  only  provide  the  relationship  formula.  The  implications  are  covered  in  [25]. 
Keywords:  Eigenvalue  problem,  numerical  approximation,  quantum  algorithms 
PACS  numbers:  03.67.Lx,  02.60.-x 
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1  Introduction 


The  study  of  the  potential  power  of  quantum  computers  has  been  a  major  theoretical  chal¬ 
lenge.  There  will  be  an  additional  incentive  to  build  a  quantum  computer  if  we  can  iden¬ 
tify  computationally  important  problems  for  which  quantum  computation  offers  significant 
speednps  over  computation  on  a  classical  computer. 

For  discrete  problems,  the  best  known  quantum  algorithms  are  due  to  Shor  and  Grover, 
see  [27,  13].  Shor’s  algorithm  for  factorization  has  an  exponential  speedup  over  all  known 
algorithms  on  a  classical  computer.  Still,  we  can  not  yet  claim  that  we  have  an  exponential 
speedup  for  this  problem,  since  the  complexity  of  factorization  on  a  classical  computer  is 
unknown.  Grover’s  algorithm  for  data  search  offers  a  quadratic  speedup. 

For  continuous  problems,  quantum  complexity  is  known  for  linear  problems  such  as  mul¬ 
tivariate  integration,  path  integration  and  multivariate  approximation,  see  [14,  15,  16,  17, 
24,  32],  For  these  problems  we  have  an  exponential  speedup  over  the  worst  case  setting,  and 
a  polynomial  speedup  over  the  randomized  setting.  The  first  quantum  study  of  a  nonlin¬ 
ear  continuous  problem  was  done  in  [20]  for  ordinary  differential  equations  with  polynomial 
speednps  over  the  classical  settings. 

The  purpose  of  this  paper  is  to  present  classical  and  quantum  complexity  results  of 
another  nonlinear  continuous  problem.  This  continuous  problem  is  quite  natural  and  com¬ 
putationally  important,  since  it  corresponds  to  the  (simplified)  univariate  Sturm-Liouville 
eigenvalue  problem.  The  Sturm-Liouville  eigenvalue  problem  is  defined  in  [10]  in  full  gener¬ 
ality.  Here  it  is  defined  as  finding  the  smallest  eigenvalue  of  the  differential  operator 

L qu{x)  :=  —u'\x)  +  q(x)  u(x)  for  x  G  (0,1), 

with  the  boundary  conditions  w(0)  =  w(l)  =  0.  We  assume  that  the  function  q  is  non¬ 
negative  and  belongs  to  the  class  C2([0,1])  of  twice  continuously  differentiable  functions 
whose  norm  ||<2j|  :=  maxi=0,i,2  max^gjoq]  \ql'l\x)\  is  bounded  by  1.  The  operator  Lg  maps 
^(M)  into  C([0,1]). 

The  Sturm-Liouville  eigenvalue  problem  has  been  extensively  studied  in  the  literature. 
The  properties  of  the  eigenvalues  and  the  eigenfunctions  are  well  known  and  so  are  numerical 
algorithms  for  approximating  them  on  a  classical  computer,  see,  e.g.  [2,  9,  10,  19,  28]. 
Nevertheless,  the  complexity  of  approximating  the  smallest  eigenvalue  in  the  worst  case  and 
randomized  settings,  as  well  as  in  the  quantum  setting,  has  not  yet  been  addressed. 

In  this  paper  we  study  classical  and  quantum  algorithms.  We  prove  bounds  on  the 
worst  case  and  randomized  complexities  on  a  classical  computer,  and  bounds  on  the  query 
complexity  and  on  the  qubit  complexity.  We  prove  that  the  complexity  in  the  classical 
settings  is  a  polynomial  in  e~l . 

We  study  the  quantum  setting  with  bit  queries  and  prove  polynomial  speedups  over  the 
classical  settings.  Bit  queries  correspond  to  approximate  computation  of  function  values, 
see  [14],  and  are  used  in  all  papers  dealing  with  the  quantum  study  of  continuous  problems. 

We  also  study  the  quantum  setting  with  power  queries.  Such  queries  are  formally  defined 
in  Section  5.2.  Here  we  only  mention  that  they  are  used  in  the  phase  estimation  algorithm, 
which  is  the  core  of  many  quantum  algorithms  including  Shor’s  and  Grover’s  algorithms. 
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Power  queries  are  controlled- WPj  queries  for  some  n  x  n  unitary  matrix  W  and  some  expo¬ 
nents  pj.  For  the  phase  estimation  algorithm,  we  have  p3  =  2-7-1  for  j  —  1,2, ,  m,  with  m 
of  order  log  £_1.  For  the  factoring  problem  of  a  large  integer  N,  Shor’s  algorithm  uses 
the  unitary  matrix  W  such  that  power  queries  can  be  implemented  by  at  most  O  (log3  IV) 
elementary  quantum  gates. 

For  the  Sturm- Liouville  eigenvalue  problem,  as  well  for  all  problems  studied  in  [25],  we 
use  power  queries  with  the  specific  unitary  matrix 

W  =  exp  (|  i  Mqj  with  i  =  V'—l,  (1) 

where  Mq  is  an  n  x  n  real  symmetric  tridiagonal  matrix  that  is  a  classical  approximation 
of  the  differential  operator  L q,  see  Section  3.2.  The  matrix  Mq  depends  on  the  values  of 
q(j / (n  +  1))  that  appear  on  the  diagonal  of  Mq  for  j  =  1,2,...,  n. 

Unitary  matrices  similar  to  (1)  play  a  key  role  in  quantum  mechanics.  They  give  the 
solution  of  the  Schrodinger  equation,  they  are  the  propagator  of  a  system  evolving  with 
Hamiltonian  Mq ,  and  are  important  in  quantum  simulation,  see  [22],  Zalka  [34]  deals  with 
their  implementation.  The  crucial  point  about  power  queries  is  that  we  can  use  W J  of  the 
matrix  W  given  by  (1)  as  one  quantum  query  for  some  j.  Hence,  lower  bound  results  for  bit 
queries  do  not  apply  to  power  queries. 

We  prove  that  in  the  quantum  setting  with  power  queries,  the  Sturm-Liouville  eigenvalue 
problem  requires  only  roughly  log  £_1  power  queries  with  the  matrix  W  of  (1).  As  shown 
in  [25],  many  computational  problems  can  be  reduced  to  the  solution  of  the  Sturm-Liouville 
eigenvalue  problem,  and  they  can  be  also  solved  in  polylog  number  of  power  queries.  The 
list  of  such  problems  include  Grover’s  search,  NP-complete  problems,  and  many  continuous 
problems.  This  proves  that  the  quantum  setting  with  power  queries  with  the  matrix  W  of 
(1)  is  exponentially  more  powerful  than  the  quantum  setting  with  bit  queries. 

We  stress  that,  contrary  to  Shor’s  algorithm,  we  do  not  know  if  power  queries  with  the 
nx  n  matrix  W  of  (1)  can  be  implemented  by  a  number  of  existing  elementary  quantum  gates 
that  is  polylog  in  n.  We  asked  a  number  of  colleagues  and  most  of  them  doubt  whether  this 
can  be  achieved.  If  this  is  indeed  the  case,  then  the  positive  results  on  the  polylog  number  of 
such  power  queries  will  be  of  only  theoretical  interest.  Still,  if  a  future  quantum  computer  is 
able  to  perform  such  power  queries  in  a  polylog  number  of,  perhaps,  more  general  elementary 
quantum  gates  or  by  some  other  quantum  devices,  the  polylog  number  of  power  queries  will 
lead  to  efficient  quantum  algorithms,  and  will  allow  us  to  solve  many  computational  problems 
exponentially  faster  than  on  a  classical  computer.  From  this  point  of  view,  we  may  interpret 
the  positive  results  on  the  number  of  power  queries  with  the  matrix  W  of  (1)  as  the  indication 
that  building  a  quantum  computer  with  such  queries  would  be  a  very  desirable  task  which 
would  give  us  a  very  powerful  computational  device. 

2  Survey  of  the  Results 

In  this  section  we  explain  our  results  in  more  technical  terms.  For  a  classical  computer,  we 
study  the  worst  case  and  randomized  settings  in  the  real  number  model  of  computation  with 
oracles,  see  [23,  30,  31].  That  is,  we  assume  that  arithmetic  operations  (addition,  subtraction, 
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multiplication,  division,  and  evaluation  of  elementary  functions),  as  well  as  comparisons  of 
real  numbers,  are  performed  exactly  with  cost  taken  as  unity.  We  also  assume  that  the 
information  about  functions  q  is  given  by  sampling  q  at  finitely  many  points  with  the  cost 
of  one  function  evaluation  taken  as  c.  Typically  c>  1. 

We  want  to  approximate  the  smallest  eigenvalue  A (q)  of  the  operator  Lg  to  within  e. 
Let  n(e)  be  the  smallest  number  of  function  values  of  q  needed  to  compute  such  an  e- 
approximation  in  a  given  setting.  The  number  n(e)  is  called  the  information  complexity.  The 
complexity ,  comp(e),  is  defined  as  the  minimal  total  cost  of  computing  an  e- approximation 
in  a  given  setting.  Obviously  we  have 

c  n{e)  <  comp(e). 

We  prove  that  in  both  classical  settings,  the  complexity  of  the  Sturm-Liouville  eigenvalue 
problem  is  polynomial  in  £~1,  or  equivalently  is  exponential  in  the  number  [log  £~1\  of 
correct  bits  of  a  computed  approximation.  More  precisely,  there  exist  positive  numbers  ctj 
independent  of  £  such  that: 


•  in  the  worst  case  setting, 


a  i  £  < 

n{£) 

< 

a2£  1/2 , 

OilC£~lJ2  < 

cornp(c) 

< 

a2c£~1/2 

+  a3  £  1/2  log  £  1 

•  in  the  randomized  setting, 

a4£-2/5  < 

n{£) 

< 

a5£~2/\ 

a±c£~2^  < 

comp(e) 

< 

a5c£~2/5 

+  ae  £~1^2  log 

The  lower  bounds  on  n(e),  and  consequently  on  cornp(c),  are  obtained  by  relating  the 
eigenvalue  problem  to  the  integration  problem  for  functions  from  the  unit  ball  of  C2([0, 1]). 
It  is  well  known  that  the  minimal  number  of  function  values  for  this  integration  problem  is 
bounded  from  below  by  roughly  e_1//2  in  the  worst  case  setting  and  by  e-2/5  in  the  randomized 
setting;  see,  e.g.,  [23,  30]  and  the  survey  of  these  results  in  [31]. 

The  upper  bounds  on  n{e)  and  comp(e)  in  the  worst  case  setting  are  obtained  by  the 
cost  of  the  classical  algorithm  that  computes  an  ^-approximation  by  the  bisection  algorithm 
of  the  Sturm  sequence  [33,  p.  300],  see  also  [11,  Ch.  5.3.4],  applied  to  an  n  x  n  matrix  which 
is  the  classical  discretization  of  the  operator  hq  with  n  =  ©(e^1/2).  The  matrix  depends 
on  n  function  values  of  q  computed  at  equidistant  points  of  [0,1].  Since  we  need  roughly 
log  e-1  bisection  steps,  and  the  cost  of  each  step  is  proportional  to  n,  the  total  cost  is  of 
order  (c  +  log  e^1)^-1/2.  Hence,  modulo  the  logarithm  of  £ -1,  the  worst  case  complexity  is 
of  order  cr1/2. 

The  upper  bounds  on  n{e)  and  comp(e)  in  the  randomized  setting  are  obtained  by  the 
following  algorithm.  We  first  approximate  the  function  q  by  a  natural  cubic  spline  q  using 
n  deterministic  sample  points  of  q  at  equidistant  points  of  [0, 1]  with  n  =  ©(e-2/5).  The 
relationship  between  the  smallest  eigenvalue  and  integration  problems,  see  Section  3,  states 
that 

A (q)  =  A (q)  +  f  (q(x)  -  q(x))  u2q(x)  dx  +  0(n”4).  (2) 

Jo 


4 


Here  Uq  is  the  normalized  eigenfunction,  f*  Uq(x )  dx  =  1,  corresponding  to  the  smallest 
eigenvalue  A  (q). 

Since  we  have  complete  information  on  the  spline  q,  we  may  approximate  A  (q)  and  Uq 
with  arbitrarily  small  error.  For  A (q),  we  achieve  an  error  of  order  £  as  in  the  worst  case 
setting,  with  cost  proportional  to  e-1/2  log  e-1.  To  obtain  an  approximation  to  Uq,  we  apply 
one  step  of  the  inverse  power  algorithm  with  an  appropriately  chosen  initial  vector.  In  this 
way  we  obtain  a  vector,  from  which  we  compute  Uq  via  piecewise  interpolation.  The  total 
cost  of  computing  A (q)  and  Uq  is  of  order  e~x!2  log  e-1. 

We  then  approximate  the  second  term  in  (2)  using  the  Monte  Carlo  algorithm  for  the 
function  (q(x)  —  q(x))u^(x)  computed  at  n  randomized  points  with  uniform  distribution  over 
[0, 1].  This  leads  to  an  ^-approximation  in  the  randomized  setting  with  cost  bounded  from 
above  by  a  quantity  proportional  to  ce~2^  +  £~1//2  log  e^1,  where  the  first  term  bounds  the 
information  cost  and  the  second  term  bounds  the  combinatorial  cost  of  the  algorithm.  Hence, 
we  have  a  sharp  estimate  on  the  randomized  information  complexity  n(s).  The  ratio  of  the 
upper  to  lower  bounds  of  the  randomized  complexity  is  roughly  at  most  e-1/10. 

In  both  classical  settings,  algorithms  for  which  we  obtain  upper  bounds  on  complexity 
require  space  of  order  This  follows  from  the  fact  that  we  need  to  work  on  n  x  n 

tridiagonal  matrices  with  n  of  order  e~1/'2 . 

We  now  turn  to  the  quantum  setting.  Quantum  algorithms  are  described  in  Section  4. 
Here  we  only  mention  that  quantum  algorithms  work  on  2V  x  2V  unitary  matrices,  where  v 
is  the  number  of  qubits.  The  qubit  complexity  is  defined  as  the  minimal  number  of  qubits 
needed  to  solve  a  problem.  Roughly  speaking,  the  qubit  complexity  corresponds  to  the 
space  complexity  for  a  classical  computer.  For  the  foreseeable  future,  qubits  will  be  a  scarce 
resource.  That  is  why  the  qubit  complexity  is  especially  important,  and  computationally 
important  problems  with  relatively  small  qubit  complexity  are  of  special  interest. 

We  prove  that  the  qubit  complexity,  compqub(e),  of  the  Sturm-Liouvillc  eigenvalue  prob¬ 
lem  is  of  order  log  £_1,  which  is  relatively  modest.  In  this  paper  log  denotes  log2.  More 
precisely,  we  prove  that 

\  log  e~l  +  0(1)  <  cornpqub (e)  <  |  log  £-1  +  0(1). 

These  bounds  hold  regardless  of  the  kind  of  queries  used.  Clearly,  the  qubit  complexity 
yields  a  lower  bound  for  the  cost  of  any  quantum  algorithm  solving  this  problem. 

We  now  turn  to  the  quantum  setting  with  bit  queries.  We  show  that  the  bit  query 
complexity  is  ©(e-1/3).  This  result  is  obtained  by  using: 

•  equation  (2)  relating  the  Sturm-Liouville  eigenvalue  problem  to  integration, 

•  a  lower  bound  on  bit  queries  for  integration,  and 

•  a  modification  of  the  classical  randomized  algorithm  described  above  that  uses  a  quan¬ 
tum  summation  algorithm  instead  of  Monte  Carlo  to  approximate  the  weighted  integral 
in  (2). 

We  now  discuss  the  quantum  setting  with  power  queries.  In  this  setting,  the  Sturm- 
Liouvillc  eigenvalue  problem  can  be  solved  using  the  well-known  phase  estimation  algorithm 
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as  a  basic  tool,  see,  e.g.,  [22,  Section  5.2],  This  algorithm  uses  power  queries  and  the 
quantum  inverse  Fourier  transform  as  its  main  ingredients.  The  power  queries  have  the  form 
controllcd-IF2J  for  j  G  N,  i.e.,  they  use  powers  of  the  matrix  W  =  exp  (|i  Mg),  with  Mq  an 
n  x  n  real  symmetric  tridiagonal  matrix  whose  diagonal  elements  depend  on  the  values  of  q. 
The  matrix  Mq  is  a  well-known  discretization  of  the  differential  operator  hq,  and  its  size  n 
depends  on  the  necessary  accuracy.  To  obtain  an  ^-approximation  we  use  n  of  order  e-1/2. 

The  phase  estimation  algorithm  uses  the  exact  eigenvector  of  Mq,  equivalently  of  W,  as 
part  of  its  initial  state,  see  [22,  Section  5.2],  Abrams  and  Lloyd  [1]  analyzed  the  case  when 
the  exact  eigenvector  is  replaced  by  an  approximate  eigenvector  and  concluded  that  as  long 
as  the  approximation  is  good  enough,  the  phase  estimation  algorithm  will  still  supply  a  good 
approximation  to  the  corresponding  eigenvalue.  Jaksch  and  Papageorgiou  [18]  proposed  an 
efficient  construction  of  an  approximate  eigenvector.  Their  idea  was  to  solve  the  problem 
with  low  accuracy  on  a  classical  computer  and  obtain  a  “short”  vector  which  approximates 
the  eigenfunction  uq  at  few  points.  Then  the  amplitudes  of  this  short  vector  are  replicated 
on  a  quantum  computer  by  the  Hadamard  transform,  which  yields  a  “long”  (vector)  state 
that  can  be  used  as  the  approximate  initial  state  in  the  phase  estimation  algorithm. 

We  show  how  the  construction  of  Jaksch  and  Papageorgiou  can  be  used  for  the  Sturm- 
Liouvillc  eigenvalue  problem.  In  this  way,  we  compute  an  e- approximation  of  the  smallest 
eigenvalue  with  probability  |  by  the  phase  estimation  algorithm  using  log  £_1  +0(1)  power 
queries.  The  algorithm  requires  an  additional  number  of  quantum  operations  at  most  of 
order  log2£_1.  This  additional  cost  is  for  the  quantum  inverse  Fourier  transform.  Finally, 
the  number  of  qubits  is  |  log  e^1  +  0(1).  A  lower  bound  on  the  number  of  power  queries 
of  order  log  £_1  has  been  proven  in  [5]. 

Comparing  these  quantum  estimates  to  the  classical  complexity  bounds  in  the  worst 
case  and  randomized  setting,  we  see  that  the  quantum  setting  with  power  queries  yields  an 
exponential  speedup  between  the  number  of  power  queries  and  the  number  of  function  values 
needed  for  the  Sturm-Liouville  eigenvalue  problem. 

Finally,  we  point  out  important  consequences  of  our  results,  which  we  study  in  detail  in 
[25].  Knowing  that  the  Sturm-Liouville  eigenvalue  problem  can  be  solved  with  polylog  power 
queries,  it  is  natural  to  study  which  computational  problems  can  be  reduced  to  this  problem. 
In  this  respect,  we  think  that  the  most  important  result  of  this  paper  is  the  formula  that 
relates  this  eigenvalue  problem  to  integration.  In  a  particular  case,  this  formula,  see  (10), 
states  that 

A(g)  =  vr2  +  1  +  2  /  (q(x)  -  |)  sin2(7nr)cfe  +  O  (||g  -  §||^)  .  (3) 

Jo 

Hence,  the  problem  of  computing  the  smallest  eigenvalue  is  equivalent,  modulo  the  second 
order  term,  to  the  weighted  integration  problem.  Since  A  (q)  can  be  approximated  with  poly¬ 
log  power  queries,  so  can  the  weighted  integral  of  q.  It  turns  out  that  many  computational 
problems  can  be  formulated  as  an  integration  problem.  Examples  include  important  discrete 
problems  such  as  Grover’s  search,  the  approximation  of  the  Boolean  mean,  and  NP-complete 
problems.  The  approximation  of  the  Boolean  mean  is  used  as  the  primary  tool  to  compute 
multivariate  integrals  and  path  integrals.  Hence,  all  these  problems  can  be  solved  by  reduc¬ 
ing  them  to  the  Sturm-Liouville  eigenvalue  problem  with  a  polylog  number  of  power  queries 
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in  the  quantum  setting.  It  is  well-known  that  Grover’s  search  and  the  approximation  of  the 
Boolean  mean  require  a  number  of  bit  queries  polynomial  in  the  problem  size,  which  in  our 
case  is  a  polynomial  in  e^1.  This  shows  that  power  queries  are  exponentially  more  powerful 
than  bit  queries,  see  [25]  for  details. 


3  Problem  Definition 


We  deal  with  functions  from  the  class 


Q  =  {  Q  :  [0, 1]  -»•  [0, 1]  |  qe  C2([0, 1])  and  ||q||  :=  max  max  |qW(x)|  <  1  }. 

For  a  function  q  G  Q,  we  consider  the  Sturm-Liouvillc  eigenvalue  problem  L qu  =  Xu  for  a 
non-zero  u,  or  equivalently 

u"{x)  —  q(x)u(x)  +  Xu(x)  =  0,  for  x  G  (0,1),  (4) 


with  the  boundary  conditions 

■u(O)  =  u(l)  =  0.  (5) 

Let  A  =  A  (q)  be  the  smallest  eigenvalue  of  (4),  (5).  Multiplying  (4)  by  u  and  integrating  by 
parts,  see  [2,  10,  28],  we  conclude  that  the  smallest  eigenvalue  satisfies 


■%) 


h1  \{u'{x))2  +  q(x)u2(x)\  dx 

min  — - - - - - - 

o f  u2(x)dx 


(6) 


where  Hq  is  the  Sobolev  space  of  absolutely  continuous1  functions  for  which  v!  G  L2([0, 1]) 
and  w(0)  =  w(l)  =  0. 

Let  uq  be  a  normalized  real  eigenfunction  corresponding  to  the  smallest  eigenvalue.  It  is 
known  that  the  eigenvalues  of  Lg  are  simple,  and  the  eigenspace  corresponding  to  A (q)  is  of 
dimension  one.  Therefore  uq  is  uniquely  defined  up  to  the  sign.  In  particular,  u2  is  uniquely 
defined.  Then  (6)  states  that 


(7) 


Observe  that  q  G  Q  implies  that  uq  G  C4([0,1]).  Since  ||g||  <  1,  and  IIuqIIlj  =  1  with 
ug(0)  =  uq(  1)  =  0,  then  \uq\x)\  are  uniformly  bounded  for  all  i  G  [0,4],  x  G  [0, 1]  and  q  G  Q, 
see  e.g.,  [10,  p.  337]. 

The  smallest  eigenvalue  A (q)  is  a  non- decreasing  function  of  q ,  i.e.,  qi(x)  <  for 

x  G  [0, 1]  implies  A(qi)  <  A (^2)-  It  is  known  that  for  q  =  c  we  have 

A(c)  =  7T2  +  c  and  uc(x)  =  a/2  sin (tix). 

1A  function  /  is  absolutely  continuous  if  and  only  if  it  can  be  written  as  f(x)  =  /( 0)  +  //  f'(t)dt  for  all 
x  €  [0, 1]. 
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This  implies  that  for  q  G  Q,  we  have  A (q)  G  [A(0),  A(l)]  =  [7 r2, 7r2  +  1], 

We  will  need  estimates  of  the  smallest  eigenvalues  and  their  eigenfunctions  for  perturbed 
functions  q.  This  is  a  classical  problem  and  many  such  estimates  can  be  found  in  the  litera¬ 
ture,  not  only  for  the  simplified  Sturm-Liouville  problem  that  we  consider  in  this  paper  but 
also  for  more  general  eigenvalue  problems.  In  our  case,  the  problem  of  perturbed  eigenvalues 
and  eigenvectors  is  well-conditioned,  since  the  differential  operator  hq  is  symmetric  and  the 
eigenvalues  of  hq  are  well  separated.  Combining  results  from  [10,  19,  29]  one  can  obtain  the 
following  estimates  for  q,  q  G  Q: 


\m-m\ 

II  'Uq  '^qlloo 

Kq) 


<  h~q\\oo  ■=  max  |g(x)  -  q(x)\, 

xe[o,i] 

=  O  (\\q  —  ?||oo)  * 

=  Ho)  +  [  (q(x)  ~q(x))u2q(x)dx  +  O  (||g-g||y  • 
Jo 


(8) 

(9) 

(10) 


We  stress  that  the  factors  in  the  big-0  notation  are  independent  of  q  and  q. 

These  relations  follow  by  elementary  arguments.  Indeed,  (8)  follows  from  (6)  by  taking 
u  =  uq,  which  leads  to  \(q)  —  X(q)  <  ||g  —  q\\oo.  By  replacing  the  roles  of  q  and  q  we  get 
\(q)  —  X(q)  <  \\q  —  q || 00,  which  implies  (8).  The  next  relation  (9)  can  be  also  proved  by 
a  matrix  approximation  to  the  operator  Ly,  which  will  be  done  in  Section  4.  Finally,  (10) 
follows  by  again  taking  u  =  uq  in  (6),  which  leads  to 


A (q)  <  X(q)  +  /  (q(x)  -  q(x))  u2q{x)  dx 


=  A (g)  +  /  (q(x)  -  q(x))  u2q{x)  dx  +  /  (q(x)  -  q(x))  [u2q{x)  -  uq(x))  dx. 


By  (9),  the  last  term  is  of  order  ||g  —  q,||^0.  Taking  u  =  uq  in  the  expression  (6)  defining  A (g), 
we  obtain 

A (q)  <  A (q)  +  /  (q(x)  -  q(x))  u2(x)  dx. 

Jo 

The  last  two  inequalities  imply  (10).  We  shall  see  later  that  the  formula  (10)  will  be  very 
useful  in  deriving  lower  bounds  for  classical  algorithms.  Note  that  if  we  take  q=\,  then  the 
formula  (10)  becomes  (3). 


4  Classical  Algorithms 

In  this  section  we  consider  classical  algorithms,  i.e.,  algorithms  on  a  classical  (non-quantum) 
computer.  These  algorithms  can  be  either  deterministic  or  randomized.  They  use  informa¬ 
tion  about  the  functions  q  from  Q  by  computing  q{ti)  for  some  discretization  points  ti  G  [0, 1]. 
Here,  i  =  1,2,...,  nq.  for  some  nq,  and  the  points  ti  can  be  adaptively  chosen,  i.e.,  ti  can  be 
a  function 

U  =  tiitx.qitx), . . .  ^ti-^qiti-x)), 

of  the  previously  computed  function  values  and  points  for  i  >  2.  The  number  nq  can  also  be 
adaptively  chosen,  see,  e.g.,  [30]  for  details. 


A  classical  deterministic  algorithm  produces  an  approximation 

o(q)  =  0(g(ti),...,g(tnJ) 

to  the  smallest  eigenvalue  A (q)  based  on  finitely  many  values  of  q  computed  at  deterministic 
points.  Let  n  =  supggQ  nq.  We  assume  that  n  <  oo.  The  worst  case  error  of  such  a 
deterministic  algorithm  0  is  given  by 

ewor(0,n)  =  sup  | A (q)  -  <f>(q)\.  (11) 

geQ 

A  classical  randomized  algorithm  produces  an  approximation  to  A (q)  based  on  finitely 
many  values  of  q  computed  at  random  points,  and  is  of  the  form 

<f>M  =  tMti  ,«)»•••,  ?(*»*,„,«)), 

where  (t>u,tijU,  and  nqjUI  are  random  variables.  We  assume  that  the  mappings 

^  1  *  ti,oo  ti{t q{t\tijd)i  •  •  •  i  qi^i—  1,oj))) 

LU  H- >  (f)U: 

LU  e- >  nq^ 

are  measurable.  Let  nq  =  E (nqtU)  be  the  expected  number  of  values  of  the  function  q  with 
respect  to  w  .  As  before,  we  assume  that  n  =  sup?eq  nq  <  oo.  The  randomized  error  of 
such  a  randomized  algorithm  0  is  given  by 

eran(0,  n)  =  sup  (E[A (q)  -  (pM)?) V2  •  (12) 

geQ 

For  simplicity  and  brevity  we  consider  the  error  of  randomized  algorithms  in  the  L2  sense. 
It  is  straightforward  to  extend  our  results  for  the  error  of  randomized  algorithms  defined  in 
the  Lp-sense  with  p  G  [1,  oo] . 

We  denote  the  minimal  number  of  function  values  needed  to  compute  an  e- approximation 
of  the  Sturm-Liouville  eigenvalue  problem  in  the  worst  case  and  randomized  settings  by 

nwor(e)  =  min {n:  3  0  such  that  ewor(0, n)  <  e}  and 
nran(e)  =  min{  n  :  3  0  such  that  eran(0,  n)  <  e  }, 


respectively. 


4.1  Lower  Bounds 

We  now  prove  lower  bounds  on  nwor(e)  and  nran(e). 

Theorem  4.1. 


nwor(e)  =  Q  (£~l/2)  ,  nran(e)  =  ft  ( e~2/ 5) 
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Proof.  Define 


F=  {/:  /eC2([0,l]),  max  (H/m,  H/'lloo,  ||/"||0o  )  <  1  })  ,  (13) 

and  consider  the  weighted  integration  problem 

1(f)  =  f  f(x)sm2(nx)  dx  V/eF. 

J  o 

It  is  well-known  that  any  algorithm  using  n  function  values  for  approximating  of  this  weighted 
integration  problem  has  worst  case  error  at  least  proportional  to  n~2  in  the  worst  case  setting, 
and  to  n-2'5  in  the  randomized  setting,  see  [23,  30] 2 . 

For  c  G  (0,  |],  consider  the  class 

Fc  =  Fn{feF  :  H/IU  <  c}.  (14) 

For  n~2  much  less  than  c,  the  proofs  for  the  class  F  can  be  used  to  deduce  the  same  lower 
bounds  on  algorithms  for  approximation  of  the  weighted  integration  problem  for  the  class  Fc. 
For  f  E  Fc  define  q  —  \  +  /.  Then  q  G  Q.  From  (3)  we  have 

A  (q)  =  vr2  +  I  +  2 1(f)  +  0(c2). 

For  any  algorithm  <f>  using  n  function  values  of  q  for  the  Sturm-Liouville  eigenvalue  problem, 
define  the  algorithm  if(f)  =  \(<f(q)  —  vr2  —  |)  for  the  weighted  integration  problem.  Then  if 
uses  n  function  values  of  /,  and 

A(g)  -  <t>(q )  =  2  (1(f)  -  if(f))  +  0(c2).  (15) 

Let  c  =  n-3/2.  Then  n~2  =  o(c),  and  therefore  the  error  of  0  is  lower  bounded  by  f }(n~2)  in 
the  worst  case  setting,  and  by  fl(n-2'5)  in  the  randomized  setting.  Hence,  the  error  of  <f  is 
at  most  e  when  n  =  Vt(e~1^2)  in  the  worst  case  setting,  and  n  =  Vt(e~2F)  in  the  randomized 
setting.  Since  this  holds  for  an  arbitrary  algorithm  0,  the  proof  is  complete.  □ 

4.2  Upper  Bounds  in  the  Worst  Case  Setting 

We  now  discuss  upper  bounds  on  nwor(e),  as  well  as  bounds  on  the  complexity  in  the  worst 
case  setting.  The  worst  case  cost  of  an  algorithm  0  using  n  function  values  is  defined  as 

costwor(0)  =  sup  (c  nq  +  mq) , 

<?gq 

where  rnq  is  the  number  of  arithmetic  operations  used  by  the  algorithm  for  a  function  q 
from  Q.  The  worst  case  complexity  compwor(e)  is  defined  as  the  minimal  cost  of  an  algorithm 
whose  worst  case  error  is  at  most  s, 

compwor(e)  =  min {  cost wor(0)  :  0  such  that  ewor(0, n)  <  e}  . 

2Formally,  these  results  are  proved  for  1(f)  =  f  )  f(x)  dx.  However,  the  same  proofs  can  be  applied  for 
the  integration  problem  with  the  weight  sin2(7ra;)  and  the  same  lower  bounds  hold. 
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Obviously,  compwor(£)  >  cnwor(£). 

We  now  discuss  the  classical  algorithm  for  the  Sturm-Liouvillc  eigenvalue  problem,  see 
e.g.,  [11,  19],  and  show  that  it  is  almost  optimal  in  the  worst  case  setting.  This  algorithm 
uses  n  =  0(£-1/2)  function  values  of  q  at  the  equidistant  points  i / (n  +  1)  for  i  —  1,  2, . . . ,  n. 
Then  the  operator  Lg  is  approximated  by  the  tridiagonal  n  x  n  matrix  Mq  of  the  form 


2  -1 

r  9(^t)  l 

-1  2  -1 

9(A) 

Mq  =  (n  +  l)2 

-1  2  -1 

+ 

9(S) 

-1  2 

L  5(A)  J 

Clearly,  Mq  is  a  symmetric  and  positive  definite  matrix.  Let  \j  =  A j(Mq)  and  Zj  =  Zj(Mq ) 
be  the  eigenvalues  and  eigenvectors  of  Mq,  i.e.,  MqZj  =  XjZj  with 


Ai  <  A2  <  •  •  •  <  An, 


where  the  vectors  Zj  are  orthogonal  and  normalized  such  that 


1 2 
\L2 


n 


n 

E-2 

k= 1 


i.fc 


1 


with  Zjtk  being  the  kth  component  of  Zj.  Note  that  we  use  the  subscript  L2  in  the  norm  of  a 
vector  to  stress  similarity  to  the  L2  norm  of  functions,  and  to  distinguish  from  the  Euclidean 
second  norm.  Clearly,  ||%||l2  =  ^||z||2- 
For  q  =  c,  it  is  known,  see,  e.g.,  [11],  that 


A  j(Mc)  =  c  +  4  (n  +  l)2si 


2  •  2 

sm 


J7T 


2  (n  +  1)7’ 


and  Zj(Mc)  =  [zj,i(Mc),  zjj2(Mc), . . . ,  Zj,n(Mc)]T  with 


Zj,k(Mc) 


2  n 

n  +  1 


1/2 

sin 


f  jkn  A 
\n+  lj  ' 


It  is  known,  see,  e.g.,  [19],  that  the  smallest  eigenvalue  Ai (Mq)  of  the  matrix  Mq  approx¬ 
imates  the  smallest  eigenvalue  A (q)  of  the  operator  hq  with  error  of  order  n~2,  i.e., 

A  (q)  -  A  .(Mg)  =  O  (n~2)  =  0(e). 

Hence,  it  is  enough  to  approximate  X\(Mq )  with  error  of  order  e.  This  can  be  achieved  by 
using  roughly  log  e~x  bisection  steps.  Each  step  consists  of  computing  the  n  terms  of  the 
Sturm  sequence,  and  this  can  be  done  in  cost  proportional  to  n.  The  total  cost  is  of  order 
(c  +  log  e~l)e~1/2.  For  details,  see  [11,  33].  Theorem  4.1  and  the  cost  of  this  algorithm  lead 
to  the  following  bounds  for  the  minimal  number  of  function  values  and  for  the  worst  case 
complexity. 
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Theorem  4.2. 


nwor(e)  =  @(e-1/2),  n(ce-1/2)  =  compwor(e)  =  0{ ce~1/2  +  W1/2log  e'1). 


Remark  4.1.  We  now  show  how  (9)  can  be  proven,  based  on  the  properties  of  the 
matrix  Mq.  First  observe  that  for  q  =  0,  the  eigenvalues  Xj(M0)  are  well  separated,  since 


Ai+i(M0)  —  \j(M0)  —  4 (n  +  l)2  sin 

>  4 (n  +  l)2  sin 


7 r 


sm 


(2  j  +  1)tt 
2(n+l)  ~“2(n+l) 

37T  7T 

sin  ■ 


2  (n  +  1)  2  (n  +  1) 


3vr2. 


For  q  E  Q,  the  Hermitian  matrix  Mq  differs  from  M0  by  the  diagonal  matrix  diagg(i/(n  +  l)) 
whose  elements  satisfy  q(i/(n+  1))  G  [0,  ||g||oo]  with  ||g||oo  <  1.  Using  the  known  estimates 
on  the  perturbed  eigenvalues  of  Hermitian  matrices,  see  [33] ,  we  have 


min 


W{Mq)  —  Aj(M0)| 


<  Iklloo 


for  all  j  =  1,2,...,  n.  Since  the  intervals  [A i(M0)  —  1,  A i(M0)  +  1]  are  disjoint,  we  conclude 
that 

|Aj(fU?)  —  Xj(M0)\  <  UglU  <  1, 

and  that 

Aj+i(Mq)  —  Xj(Mq)  >  Aj+i(M0)  —  Xj(M0)  —  2  «  37r2  —  2. 


Define 


W<7,n  — 


Un 


n  +  1 


>  uq 


n 


n+1 


where  uq  is  the  normalized  real  eigenfunction  corresponding  to  the  smallest  eigenvalue.  Then 
II || Z/2  =  1  +  o(l).  We  normalize  uq^n  and  obtain 


1 


Hq,n  — 


■  U, 


\uq,n\\L2 


q,n • 


As  mentioned  in  Section  3,  the  eigenfunction  uq  is  defined  uniquely  up  to  its  sign.  Ob¬ 
viously,  the  same  is  true  for  the  eigenvector  Zi(Mq).  We  choose  the  signs  of  uq  and  Zi(Mq ) 
such  that 

\\Uqtn  -  Zi(Mq)\\L2  <  \\uq,n  +  Zi  (Mq)  ||  l2  ■ 

All  the  components  of  the  vector 


T]n  ■ —  MqUq^n  ^{q)^q,n 

are  of  order  n-2,  and  therefore  1 1 77n 1 1 =  0(n~2).  From  the  a  posteriori  error  estimate,  see 
[33,  p.  173],  we  conclude  that 

| \uq>n  -  Zi(Mq)\\L2  =  0(n~2)  Wq  G  Q 
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with  the  factor  in  the  big-0  notation  independent  of  q.  Note  also  that 

Mqtlq  n  MqUq  n  \(yq)Uqn  T  rn, 

with  ||rn||ia  =  0(||g  -  g||oo).  Hence 

\\uq,n  ~  Zi(Mq)\\L2  =  0(\\q  —  gHoo  +  n~2). 

Finally,  we  have 

\\Uq,n  -Mff.nllLa  =  IK,n  ~  ^1  (Mq)  +  ZX  (Mq)  -  Uq>n  ||La  =  O  (n'2  +  ||  q  -  q\\  <*,) . 
Letting  n  tend  to  infinity,  we  conclude  that 

Wuq-uq\\L2  =  0(\\q-q\U. 

Since  both  uq  and  Ug  satisfy  (4)  for  (q,  A (q))  and  (q,  A (q)),  respectively,  we  have 

Uq(x)  -  uq(x)  =  (q(x)  -  A (q))(ug(x)  -  Ug(x))  +  Ug(x)  (( q(x )  -  q(x))  -  (A (q)  -  A (q)))  . 
Therefore 

IK^X-11^2  =  0(\\q  —  g||oo). 

This  and  the  fact  that  u  —  Ug  vanishes  at  0  and  1  imply 

W'U'q  ^g||oo  0(||^  1 1  oo )  ? 


as  claimed.  □ 

4.3  Upper  Bounds  in  the  Randomized  Setting 

We  now  turn  to  the  randomized  setting.  The  cost  of  a  randomized  algorithm  0,  using 
n  =  supygq  E (nqtU1)  <  oo  randomized  function  values,  is  now  defined  as 

costran(0)  =  sup  (E  (c  nqiUJ  +  mq^)2)1/2  , 
qe  Q 

where  rnqijj  is  the  number  of  arithmetic  operations  used  by  the  algorithm  for  a  function  q 
from  Q  and  a  random  variable  u.  The  randomized  complexity 

compran(e)  =  min  {  costran(0)  :  (j)  such  that  eran(0,  n)  <  e}  , 

is  the  minimal  cost  of  an  algorithm  whose  randomized  error  is  at  most  e.  Obviously, 
cornpran(c)  >  cnran(e). 

We  now  derive  upper  bounds  on  nrail(e)  and  compran(e)  by  presenting  a  randomized  algo¬ 
rithm  that  depends  on  a  number  of  parameters.  Then  we  find  the  values  of  these  parameters 
for  which  the  randomized  error  is  e.  We  first  compute  m  + 1  function  values  of  q  at  determin¬ 
istic  points  i/m,  for  i  =  0, 1, . . .  ,m,  and  construct  a  cubic  natural  spline  gcub  interpolating 
q  at  these  points,  see  e.g.,  [7]  for  information  about  cubic  splines.  It  is  well  known  that  this 
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can  be  done  with  cost  proportional  to  m,  and  ||g  — gCub||oo  =  0(m~2).  The  function  gcub  does 
not  have  to  be  non- negative.  Since  q  >  0  then  q  =  qcub  +  c  >  0  with  a  constant  c  =  0(m~2). 
We  have  q  G  Q  and  ||g  —  g||oo  =  0(m~2).  We  apply  the  formula  (10)  for  the  function  q  and 
obtain 

A (q)  —  A (q)  =  /  (q(x)  —  q(x))  Uq(x)  dx  +  O  (m-4)  .  (16) 

Jo 

This  suggests  that  we  can  improve  the  accuracy  of  approximating  A (q)  —  A (q)  by  using  the 
classical  Monte  Carlo  algorithm  applied  to  the  first  term  of  the  right  hand  side  of  (16). 
We  will  need  to  know,  at  least  approximately,  the  eigenvalue  A (q)  and  the  eigenvector  Uq. 
Suppose  we  approximate  A  (q)  by  \q  with  the  worst  case  error 

sup  |A(g)  -  Xq\  <  5i,  (17) 

<je  Q 

and  the  eigenfunction  Uq  by  Zq  with  the  worst  case  error 

sup  || Uq  -  Zq\\L2  <  S2.  (18) 

<?gq 

Assume  for  a  moment  that  A q  and  Zq  have  been  computed.  For  a  function  v,  define 
fv(x)  =  (g(x)  -  q(x))v2(x)  and  I(fv)  =  fv(x)  dx. 

The  randomized  algorithm  (j)  based  on  the  Monte  Carlo  with  k  randomized  samples  takes 
the  form 

<t>M)  =  X<i  +  (^(^>)  _  fe>)  j-|(^>), 

K  i= i  ^  ' 

where  XjtUJ  are  independent  and  uniformly  distributed  numbers  from  [0, 1].  Here  u  represents 
a  random  element.  We  have 


ia(?)  -  mq)\  <  \m  -  a,i  + 1 i{u  -  Hu i 


Zq  \Xj,U!j 


3= 1 


+  0{m 


-4n 


Clearly, 


Hfnq)  -I(fzq)  \  <  /  I q(x)  -  q{x)  \  \ u^x)  -z2q(x)\  dx  =  0{m  2  S2). 


Since  ||/2Jl2  =  0(rri  2),  the  well  known  formula  for  the  randomized  error  of  Monte  Carlo 
yields  that 


1/2 


Zq  \Xj,u) 


3= 1 


W2Zq)-nh9)Y/2 

k1/ 2 


=  O  (m-2*;-1/2)  . 


We  have  obtained  the  bound 

eran(0,  n)  =  O  (<5i  +  m_2<52  +  m~2 A;-1/2  +  m-4) 
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on  the  randomized  error  of  0.  Hence,  to  guarantee  error  at  most  e,  it  is  enough  to  take 


hi  =  ©(e),  m  —  k  —  0(e  2^5)  and  h2  =  ©(e1^5)- 


We  now  explain  how  to  achieve  (17)  and  (18).  To  get  Xq  approximating  A (q)  with  error 
of  order  e,  we  approximate  the  operator  L q  by  the  matrix  Mq  as  in  the  worst  case  setting, 
now  with  n  =  Q(e~1/2).  Then  A (q)  —  Ai (Mq)  =  0(n~2)  =  0(e),  and  we  compute  A q  as  an 
e- approximation  of  Ai  (Mq)  as  for  the  worst  case  setting.  This  can  be  done  with  cost  of  order 
e-1/2  function  values  of  q,  and  of  order  e_1,/2  log  e_1  arithmetic  operations.  Since  the  cost 
of  computing  one  function  value  of  q  is  of  order  1,  the  total  cost  of  computing  A q  is  of  order 
e-1/2  log  e-1. 

To  get  Zq  approximating  Uq  with  error  of  order  e1/5  we  proceed  as  follows.  Consider  the 
eigenvector  z\(Mq)  of  the  matrix  Mq ,  with  n  not  yet  specihed.  By  Remark  4.1,  we  have 

\\uq,n  -  Zi(Mq)\\L2  =  0(n~2).  (19) 


We  approximate  the  smallest  eigenvalue  Ai  (Mq)  by  A,  with  error  5.  This  can  be  achieved 
with  cost  of  order  n  log  W1.  Without  loss  of  generality  we  assume  that  A  ^  Ai  (Mq).  Indeed, 
we  can  check  this  condition  by  computing  det (Mq  —  A I)  and  if  this  determinant  is  zero  we 
perturb  A  a  little.  Then  the  matrix 


A  =  (Mq  -  XI) 


is  non-singular  and  its  eigenvalues  are  (i3  =  (Xj(Mq)  —  A)”1.  Note  that  \(5\\  >  and 
f3j  =  0(1)  for  j  >  2.  For  the  j th  vector  e3  —  [0, . . . ,  0, 1,  0, . . . ,  0]T  with  1  in  the  j'tli  position, 
dehne 

Xj  =  Aej. 


We  can  compute  xj  with  cost  of  order  n  by  solving  the  tridiagonal  linear  system  (Mq-XI)xj  = 
e3.  Then  we  compute 


\x 


toll  2 


=  max 

j= 1,2, ...,n 


\X 


ill  2, 


and 

'y  1 1  of*  'll  Of*  • 

^  \\,l3o\\2  ^3  o* 

Observe  that  the  cost  of  computing  z  is  of  order  n2. 

Since  {n_1/2^(Mg)}"=1  is  orthonormal,  we  have  ||a:j||2  =  Pe(eji  zt(Mq))2  and 
||n-1/2  Zi(Mq)\\l  =  1  =  Y3'j=\{e3i  n~1^2  zi(Mq))2.  Hence,  there  exists  an  index  j  such  that 


(e.j,n  1/2  Z\ (Mq))2  >  n  ', 


and  therefore 

We  have 


\x 


J0 


2  > 


\xj\\2  >  5  1n  1//2. 


(Mq  -  X1(Mq)I)z  =  (Mq  -  A I)z  +  (A  -  A x(Mq))z 


- — 7—e3o  +  (A  —  Ai  (Mq))z, 

Xjo  I  2 
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and  therefore 


|| (Mg  -  A  1(M*)i>||2  <  6y/E  +  6. 

From  [33,  p.  173],  we  conclude  that  ||n-1//2  z\(Mq)  —  z||2  =  0(8  y/ n ),  and 

\\zi(Mq)  —  \fnz\\L2  —  0(6y/n).  (20) 

We  are  finally  ready  to  define  Zq  by  piecewise  linear  interpolation  from  the  successive 
components  of  the  vector  yfnz  =  [zi,  z2, . . . ,  zn]T .  More  precisely,  for  j  =  0, 1, . . . ,  n  let 
tj  =  j/(n  +  1).  For  t  G  [tj,tj+ 1],  we  set 

zq(t)  =  Zj(l  -  (n  +  l)t  +  j)  +  zj+i((n  +  l)t  —  j) 


with  zq  =  zn+i  =  0. 

We  need  to  estimate  Uq  —  Zq  in  the  L2  norm.  Observe  that  for  t  G  [tj,tj+ 1]  we  have 
Uq(t)  =  Uq(tj)(l  ~  (n+l)t  +  j)  +  Uq(tj+1)((n  +  l)t  -  j)  +  0(n~ 2) 
since  Uq  G  Q.  Therefore 

\Uq(t)  -  Zg(t)\  <  | Uq(tj)  -  Zq(tj)\  +  \uq(tj+1)  -  Zq(tj+1)\  +  0(n~2). 

This  yields 

n  ^  ftj+ 1 

IK'-^-IlL  =  5Z  /  (Uq(t)  -  Zg(t)f  dt 

j= 0  i 

(  1  ” 

=  0  r— y  Kfo) _  %(^))2  +  n ^ 

\  3=0 

Hence, 

\\uq-Zq\\L2  =  o  (\\uqtn  -  y/nz\\L2  +  n~2) 

Since  || UqiU  —  y/nz\\L2  <  \\uq,n  —  zi(Mq)\\L2  +  \\zi(Mq)  —  y/nz\\L2,  we  use  (19)  and  (20)  to 
see  that 

\\uq-zq\\L2  =  0(5y/n  +  n~2). 

For  6  =  n-5/2  we  obtain 

\\uq-  zqh2  =  0(n~2). 

Setting  n  =  ©(e^1/10)  we  obtain  (18)  with  S2  =  ©(e1/5).  The  cost  of  computing  zq  is  of  order 

n2  =  ©(e-1/5). 

Theorem  4.1  and  the  cost  of  this  randomized  algorithm  lead  to  the  following  bounds  on 
the  minimal  number  of  function  values  and  the  randomized  complexity. 

Theorem  4.3. 

nrail(e)  =  0(e_2/5),  n(c£~2/5)  =  compran(e)  =  0(ce~2'b  +  £~x'2  log  £~1). 
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5  Quantum  Setting 

We  now  turn  our  attention  to  the  quantum  setting.  In  this  setting,  we  are  using  hybrid 
algorithms  that  are  combinations  of  classical  algorithms  using  function  values,  as  explained 
in  the  previous  sections,  and  quantum  algorithms  which  we  now  describe.  A  quantum 
algorithm  applies  a  sequence  of  unitary  transformations  to  an  initial  state,  and  the  final 
state  is  measured,  see  [3,  8,  14,  22]  for  the  details  of  the  quantum  model  of  computation. 
We  briefly  summarize  this  model  to  the  extent  necessary  for  this  paper. 

The  initial  state  j^o)  is  a  unit  vector  of  the  Hilbert  space  7d„  =  C2  (8)  •  •  •  <S>  C2,  v  times, 
for  some  appropriately  chosen  integer  v,  where  C2  is  the  two  dimensional  space  of  complex 
numbers.  Obviously,  the  dimension  of  7 ~LV  is  2".  The  number  v  denotes  the  number  of  qubits 
used  in  quantum  computation. 

The  final  state  \ij})  is  also  a  unit  vector  of  7 ~tv  and  is  obtained  from  the  initial  state  l^o) 
by  applying  a  number  of  unitary  2V  x  2V  matrices,  i.e., 

If/’)  :=  UtQyUt-iQy  ■  ■  ■  UiQyUoI^o)-  (21) 

Here,  Uq,  U\, . . .  ,Ut  are  unitary  matrices  that  do  not  depend  on  the  input  function  q.  The 
unitary  matrix  Qy  with  Y  =  . . .  ,q(tn)\  is  called  a  quantum  query  and  depends  on 

n,  with  n  <  2U,  function  evaluations  of  q  computed  at  some  non- adaptive  points  tt  e  [0,1]. 
The  quantum  query  Qy  is  the  only  source  of  information  about  q.  The  integer  T  denotes 
the  number  of  quantum  queries  we  choose  to  use. 

At  the  end  of  the  quantum  algorithm,  a  measurement  is  applied  to  its  final  state  The 
measurement  produces  one  of  M  outcomes,  where  M  <  2U .  Outcome  j  G  {0, 1, . . . ,  M  —  1} 
occurs  with  probability  py{j ),  which  depends  on  j  and  the  input  Y.  For  example,  if  M  =  2V 
and  the  final  state  is  | tjj)  =  cj\j)!§  whh  \ci\2  =  1)  then  a  measurement  in 

the  computational  orthonormal  basis  (|j)}  produces  the  outcome  j  with  probability  py  (j )  = 
\cj\2.  Knowing  the  outcome  j,  we  compute  an  approximation  A y(j)  of  the  smallest  eigenvalue 
on  a  classical  computer. 

In  principle,  quantum  algorithms  may  have  many  measurements  applied  between  se¬ 
quences  of  unitary  transformations  of  the  form  presented  above.  However,  any  algorithm 
with  many  measurements  and  a  total  of  T  quantum  queries  can  be  simulated  by  a  quantum 
algorithm  with  only  one  measurement  at  the  end,  for  details  see  e.g.,  [14]. 

We  stress  that  classical  algorithms  in  floating  or  fixed  point  arithmetic  can  also  be  written 
in  the  form  of  (21).  Indeed,  all  classical  bit  operations  can  be  simulated  by  quantum  com¬ 
putations,  see  e.g.,  [4],  Classically  computed  function  values  will  correspond  to  bit  queries 
which  we  discuss  in  Section  5.2. 

In  our  case,  we  formally  use  the  real  number  model  of  computation.  Since  the  Sturm- 
Liouvillc  eigenvalue  problem  is  well  conditioned  and  properly  normalized,  we  obtain  practi¬ 
cally  the  same  results  in  floating  or  fixed  point  arithmetic.  More  precisely,  it  is  enough  to  use 
0(log  e~l)  mantissa  bits,  and  the  cost  of  bit  operations  in  floating  or  fixed  point  arithmetic 
is  of  the  same  order  as  the  cost  in  the  real  number  model  multiplied  by  a  power  of  log  £ -1. 

Hence,  a  hybrid  algorithm  may  be  viewed  as  a  finite  sequence  of  algorithms  of  the  form 
(21).  It  is  also  known  that  if  we  use  finitely  many  algorithms  of  the  form  (21)  then  they  can 
be  written  as  one  quantum  algorithm  of  the  form  (21),  see  [14,  15]. 
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That  is  why  an  arbitrary  hybrid  algorithm  in  the  quantum  setting  is  of  the  form  (21). 
This  is  important  when  we  want  to  prove  lower  bounds  because  it  is  enough  to  work  with 
algorithms  of  the  form  (21).  For  upper  bounds,  it  seems  to  us  more  natural  to  distinguish 
between  classical  and  quantum  computations  and  charge  their  cost  differently.  The  cost 
of  classical  computations  is  defined  as  before  whereas  the  cost  of  quantum  computations  is 
defined  as  the  sum  of  the  number  of  quantum  queries  multiplied  by  the  cost  of  one  query, 
and  the  number  of  quantum  operations  besides  quantum  queries,  ft  will  be  also  important 
to  indicate  how  many  qubits  are  used  by  the  quantum  computations. 

We  now  define  the  error  in  the  quantum  setting.  In  this  setting,  we  want  to  approximate 
the  smallest  eigenvalue  A (q)  with  a  probability  p  >  \.  For  simplicity,  we  take  p  —  |  for 
the  rest  of  this  section.  As  it  is  common  for  quantum  algorithms,  we  can  achieve  an  e- 
approximation  with  probability  arbitrarily  close  to  1  by  repetition  of  the  original  quantum 
algorithm,  and  by  taking  the  median  as  the  final  approximation. 

The  local  error  of  the  quantum  algorithm  with  T  queries  that  computes  A y(j)  for  the 
function  q  G  Q  and  the  outcome  j  G  {0, 1, . . . ,  M  —  1}  is  defined  by 

e(A Y,  T)  =  min  j  a  :  ^  pY(j)  >  f  j- 

j:  |A(g)-Ay(j)l  <a 

This  can  be  equivalently  rewritten  as 

e(A y,T)  =  min  max  |A(g)  —  Ay(j)|, 

A:n(A)>\  i£A 

where  A  C  {0, 1, . . . ,  M  —  1}  and  p(A)  =  ^jeApY(j)- 

The  worst  probabilistic  error  of  a  quantum  algorithm  A  with  T  queries  for  the  Sturm- 
Liouvillc  eigenvalue  problem  is  defined  by 

equant(A,  T )  =  sup  |  e(Ay,  T) :  Y  =  [, q(ti ), . . . ,  q(tn)],  U  G  [0, 1],  for  q  G  Q  j.  (22) 

5.1  Bit  Queries 

Quantum  queries  are  important  in  the  complexity  analysis  of  quantum  algorithms.  A  quan¬ 
tum  query  corresponds  to  a  function  evaluation  in  classical  computation.  By  analogy  with 
the  complexity  analysis  of  classical  algorithms,  we  analyze  the  cost  of  quantum  algorithms  in 
terms  of  the  number  of  quantum  queries  that  are  necessary  to  compute  an  e- approximation 
with  probability  Clearly,  this  number  is  a  lower  bound  on  the  quantum  complexity,  which 
is  defined  as  the  minimal  total  cost  of  a  quantum  algorithm  that  solves  the  problem. 

Different  quantum  queries  have  been  studied  in  the  literature.  Probably  the  most  com¬ 
monly  studied  query  is  the  bit  query.  For  a  Boolean  function  /  :  {0, 1, . . . ,  2m  —  1}  — >  {0, 1}, 
the  bit  query  is  defined  by 

Qf\j)\k)  =  mk®m). 

Here  v  —  m  +  1,  |j)  G  Hm,  and  | k)  G  Hi  with  ©  denoting  the  addition  modulo  2.  For  real 
functions,  such  as  functions  q,  the  bit  query  is  constructed  by  taking  the  most  significant 
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bits  of  the  function  q  evaluated  at  some  points  tr  More  precisely,  as  in  [14],  the  bit  query 
for  q  has  the  form 

Qq\j)\k)  =  \j)\k®P(q(T(J)))), 

where  the  number  of  qubits  is  now  v  =  m'+m"  and  |  j)  G  TLm' ,  | k)  G  Hm"  with  some  functions 
(5  :  [0, 1]  — ■>  {0, 1, ... ,  2m"  —  1}  and  r  :  {0, 1, . . . ,  2m'  —  1}  — ■>  [0, 1].  Hence,  we  compute  q  at 
tj  =  r(j)  G  [0, 1]  and  then  take  the  m"  most  significant  bits  of  q(tj)  by  (3(q(tj )),  for  details 
and  a  possible  use  of  ancilla  qubits  see  again  [14]. 

Using  bit  queries,  the  well  known  quantum  algorithm  of  Grover  [13]  requires  ©(IV1/2) 
queries  for  searching  an  unordered  database  of  N  items.  Similarly,  the  quantum  summation 
algorithm  of  Brassard  et  al.  [6]  computes  the  mean  of  a  Boolean  function  defined  on  the 
set  of  N  elements  with  accuracy  £  and  probability  |  using  of  order  min{./V,  £-1}  bit  queries. 
Both  algorithms  are  optimal  modulo  multiplicative  factors  in  terms  of  the  number  of  bit 
queries. 

The  quantum  summation  algorithm  can  be  also  used  for  the  approximate  computation 
of  the  mean  of  a  real  function  /  :  [0, 1]  — >  M  with  \f(x)\  <  M  for  all  x  G  [0, 1],  see  [14,  23]. 
More  precisely,  if  we  want  to  approximate 

N—l 

SnU)  := 

1=0 

for  some  Xj  G  [0,1]  and  N,  then  the  quantum  summation  algorithm  QS N(f)  approximates 
Sjv(/)  such  that 

|Sjv(/)  -  QSJV(/)|  <  £  with  probability  |  (23) 

using  of  order  mm(N,Me~1)  bit  queries,  min(iV,  Me~l)  log  N  quantum  operations,  and 
log  N  qubits. 

Bit  queries  have  been  also  used  for  a  number  of  continuous  problems  such  as  multivariate 
and  path  integration,  multivariate  approximation,  and  ordinary  differential  equations.  Tight 
bit  query  complexity  bounds  are  known  for  a  number  of  such  problems,  see  [14,  15,  16,  17, 
20,  24,  32], 

In  particular,  Novak  [24]  proved  that  for  the  integration  problem  f(x)  dx  for  functions 
/  from  the  class  F  given  by  (13),  the  bit  query  complexity  is 

nbit— query  (£;INTf)  =  0(£-V3)_  (24) 

Here  and  elsewhere  by  the  bit  query  complexity  we  understand  the  minimal  number  of  bit 
queries  needed  to  compute  an  e-approximation  to  a  given  problem  with  probability  |.  In 
particular,  nblt‘query  (e)  denotes  the  bit  query  complexity  of  the  Sturm-Liouville  eigenvalue 
problem. 

Based  on  the  result  (24)  of  Novak  and  the  relationship  between  the  Sturm-Liouville 
eigenvalue  problem  with  integration,  we  now  prove  the  following  theorem. 

Theorem  5.1. 

nbit-query(£)  =  O^"1/3). 
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Proof.  We  first  prove  that  the  bit  query  complexity  for  the  weighted  integration  problem 
for  the  class  Fc  given  by  (14)  is  of  the  same  order  as  for  integration  for  the  class  F, 

nbit-qUery(^INTFj  =  Q(£-l/3)_  (25) 


The  upper  bound  follows  from  (24).  To  prove  the  lower  bound,  we  use  the  standard  proof 
technique  of  reducing  the  integration  problem  to  the  mean  Boolean  summation  problem  for 
which  a  lower  bound  on  bit  queries  is  known. 

Assume  then  that  we  use  an  arbitrary  quantum  algorithm  with  k  bit  queries  that  com¬ 
putes  an  ^-approximation  with  probability  |  for  the  integration  problem  over  the  class  Fc. 
Without  loss  of  generality  we  assume  that  k~2  <  c. 

Consider  the  function  h(x)  =  ax3  (l  —  a:)3  for  x  E  [0, 1]  and  h(x)  —  0  for  x  >  1.  Here,  a 
is  a  positive  number  chosen  such  that  h  E  F  with  F  given  by  (13).  For  j  —  0, 1, . . . ,  N  —  1, 
with  N  >  k,  define  hj(x)  =  N~2h(N(x  —  j/N)).  Clearly,  hj  E  F  and  the  support  of  hj  is 
(j  /N,  ( j  +  1  )/N .  Observe  that  H^Hoo  <  N~2 .  Hence  hj  E  Fc.  We  also  have  hj(x )  dx  = 
N~3  fo  h(x)  dx.  For  an  arbitrary  Boolean  function  B  :  {0, 1, . . . ,  N  —  1}  — >  {0, 1},  define  the 
function 

JV-l 

Ib(x)  =  ^2  B(J)hj(x)  y  x  E  [0,1]. 

3=0 


Then  fB  E  Fc  and 


fB(x)  dx 


foh(x)  dx  lV/h:,b 

N2  N  ^  J 

3=0 


Hence,  modulo  the  factor  of  order  N~2,  the  computation  of  the  Boolean  mean  is  reduced  to 
the  integration  problem.  Note  that  /s(t)  =  B(j)hj(t )  if  t  E  [j/N,  ( j  +  1  )/N],  and  sampling 
of  fB  is  equivalent  to  sampling  of  B.  From  [21]  we  know  that  H(A;_1)  is  a  lower  bound 
for  the  error  of  the  quantum  approximation  of  the  Boolean  mean,  with  k  bit  queries,  and 
probability  |,  where  N  >  j3k  for  some  positive  f3.  Letting  N  =  \/3k\,  we  conclude  that  the 
corresponding  lower  bound  on  the  integration  problem  over  the  class  Fc  is  Q(k~3).  Hence  to 
achieve  the  error  e  we  must  have  k  =  h2(e-1/3),  as  claimed  in  (25). 

The  same  proof  techniques  allows  us  to  consider  the  classes  Fcu\  with  varying  c(e),  even 
with  c(e)  tending  to  zero,  although  not  too  fast.  We  have 


^bit-query 


(e,INTi 


c(e) ' 


0(e~1/3) 


if  lim  c(e)  £  2//'3 
£ — >0 


(X). 


(26) 


We  now  turn  to  the  Sturm-Liouville  eigenvalue  problem.  As  in  the  proof  of  Theorem  4.1, 
for  f  E  Fc  with  c  E  (0,  |],  we  define  q  =  |  +  /  and  consider  an  arbitrary  quantum  algorithm  (p 
that  uses  k  quantum  bit  queries  and  computes  an  ^-approximation  of  the  smallest  eigenvalue 
with  probability  |.  Then  ip(f)  =  \{<j>{q)  —  vr2  —  |)  is  a  quantum  algorithm  for  approximating 
the  integration  problem  over  the  class  Fc.  We  have 


=  ||(A(g)  -</>(g))  +  0(c2)|  <  i£  +  0(c2). 

Take  now  c  =  c(e)  =  Q(e2^3~s)  with  <5  E  (0,  |).  Then 

[/(/)—  ip(f)  |  <  \  e  +  0(e:4//3_25)  =  ie(l  +  o(l))  <  £  for  small  e. 
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Hence,  the  quantum  error  of  0  with  probability  |  is  £,  and  0  uses  k  bit  queries.  Due  to  (26), 
we  have  k  =  fl(e_1//3)  which  completes  the  proof.  □ 


We  now  derive  upper  bounds  on  the  bit  query  complexity  nblt"query(e)  and  on  the  total 
quantum  complexity  compblt"quant(e).  The  total  quantum  complexity  is  defined  as  the  mini¬ 
mal  cost  of  a  hybrid  algorithm  that  solves  the  Stnrm-Lionvillc  eigenvalue  problem  with  error 
at  most  e  and  probability  |.  The  hybrid  algorithm  may  require  some  classical  computations 
and  the  use  of  function  valnes  and  the  cost  of  them  is  defined  just  as  before.  It  may  also 
require  some  quantum  computations  and  the  cost  of  them  is  defined  as  the  sum  of  the  num¬ 
ber  of  bit  queries  multiplied  by  the  cost  of  one  such  query  plus  the  number  of  additional 
quantum  operations.  The  cost  of  one  bit  query  is  denoted  by  Cbit- 

We  present  a  hybrid  algorithm,  which  will  be  a  combination  of  the  classical  algorithm 
from  Section  4  and  the  quantum  summation  algorithm  QSN  for  a  properly  chosen  N.  We 
proceed  as  in  Section  4  and  use  the  same  notation.  From  (16),  (17),  and  (18),  we  have 

A (q)  =  A q  +  f  ( q(x )  -  q(x))zq{x)  dx  +  O(<0  +  m~282  +  m~4)  (27) 

Jo 

with  Ai ,  52  and  m  to  be  specified  later.  Let 

/ (x)  =  ( q(x )  -  q(x))zg(x)  x  G  [0, 1]. 

Observe  that  f(x )  =  0(m~2),  and  f(x)  depends  on  q(x),  and  q(i/m )  for  i  =  0,1,..., m, 
which  are  used  in  the  construction  of  q.  Furthermore,  we  can  compute  f(x)  by  computing 
one  function  value  q(x)  and  one  function  value  of  the  already  computed  functions  q  and  zq 
at  x.  We  approximate  fQ  f(x)  dx  by 

j= 0  v  7 

with  N  =  (m  +  1  )k,  where  the  parameters  m  and  k  will  be  specified  later.  Since  /  is  twice 
continuously  differentiable  and  f"(x)  is  uniformly  bounded  on  the  subintervals  (i/m,  (i  + 
l)/m)  for  i  =  0, 1, . . . ,  m  —  1,  it  is  easy  to  see  that 

J  !(x)dx  -  s „(/)  = 

We  dehne  N  such  that  N~2  is  of  order  £. 

We  now  apply  QS N(f)  algorithm  to  compute  an  ©(^-approximation  with  probability  | 
t°  SN(f),  or?  equivalently  to  f(x)  dx.  To  do  it,  we  need  to  use  the  bit  query  Q f  for  the 
function  /,  although  so  far  we  assumed  that  we  can  use  only  bit  queries  Qq  for  the  functions  q 
from  Q.  This  problem  is  resolved  in  Section  2  of  [15]  where  it  is  shown  that  algorithms  using 
the  bit  query  Qf  can  be  simulated  by  algorithms  using  bit  queries  Qq  at  the  expense  of 
multiplying  the  number  of  bit  queries  by  a  factor  of  2. 

From  this  and  (23)  with  M  =  0(m-2),  we  conclude  that  its  is  enough  to  perform  of  order 
min(£-1/2,  m_2£-1)  bit  queries,  min(£-1/2,  m_2£_1)  log  e~l  quantum  operations,  and  using 
of  order  log  £-1/2  qubits. 
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We  finally  approximate  A (q)  by  the  following  algorithm 

0(g)  =  A,-  +  QS  N(f).  (28) 

This  algorithm  differs  from  the  randomized  algorithm  of  Section  4  since  we  now  apply  the 
QS^y  quantum  algorithm  instead  of  Monte  Carlo  to  approximate  fQ  f(x)  dx.  Its  error  is 
clearly  of  the  form 

gbit-quant^  T)  =  Q  +  m~^2  +  m“4  +  e)  . 

To  guarantee  that  this  error  is  at  most  e,  we  take 

<5i  =  ©(e),  m  =  0(e_1A3),  k  =  ©(e^1/6)  and  <f2  =  ©(e1^3). 

Using  the  cost  analysis  of  Section  4  and  the  results  of  this  section,  we  conclude  the  following 
theorem. 

Theorem  5.2.  The  Sturm- Liouville  eigenvalue  problem  can  be  solved  in  the  quantum  setting 
with  bit  queries  by  the  algorithm  0  defined  by  (28).  This  algorithm  approximates  the  smallest 
eigenvalue  \(q)  with  error  at  most  £  and  probability  |  using  of  order 

•  e ^1/3  bit  queries  and  function  values, 

•  e^1/3  log  e-1  quantum  operations, 

•  £~li2  log  £~l  classical  operations, 

•  log  e-1  qubits. 

Furthermore, 

nbit-query  = 

and 

U(cbite_1/3)  =  compblt_query(e)  =  O  (( c  +  cbit)e'1/3  +  e_1/2  log  e_1)  . 

Hence,  we  have  a  sharp  bound  of  order  £-1/3  on  the  number  of  bit  queries  whereas  the 
upper  bound  on  the  total  cost  depends,  as  in  the  worst  case  and  randomized  settings,  on 
£  1/'^  log  £_1,  which  is  the  cost  of  classical  computations. 

5.2  Power  Queries 

In  this  subsection  we  study  power  queries.  We  formally  define  them  as  follows.  For  some 
problems,  a  quantum  algorithm  can  be  written  in  the  form 

10)  :=  ■  ■  ■  U, Wx U0 \ 0O ) ■  (29) 

Here  U{, ... ,  Urn  denote  unitary^ matrices  independent  of  the  function  q,  just  as  before, 
whereas  the  unitary  matrices  Wj  are  of  the  form  controlled- Wj,  see  [22,  p.  178].  That 
is,  Wj  =  WPj  for  an  n  x  n  unitary  matrix  W  that  depends  on  the  input  of  the  computational 
problem,  and  for  some  non- negative  integers  pj,  j  =  1,2, . . . ,  rri.  Without  loss  of  generality 
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we  assume  that  n  is  a  power  of  two.  Let  {\yk)}  be  orthonormalized  eigenvectors  of  W, 
W \yk)  =  otk\Hk)  with  the  corresponding  eigenvalue  ak,  where  \ay\  =  1  and  ak  =  elAfc  with 
A k  G  [0,  2n)  for  k  —  1,  2, . . . ,  n.  For  the  unit  vectors  \xe)  =  at |0)  +  |1)  G  C2,  £  —  1,  2, . . . ,  r, 

the  quantum  query  Wj  is  defined  as 

Wj  |»i)|x2)  •  •  •  \xr)\yk)  =  \x1)\---\xj_l)(aj\0)+pje1'yPjXk\l)^\x:j+l)---\xr}\yk).  (30) 


Hence,  Wj  is  a  2V  x  unitary  matrix  with  v  =  r  +  log  n.  We  stress  that  the  exponent  pj 
only  affects  jthe  power  of  the  complex  number  e17Afc. 

We  call  Wj  a  power  query  since  they  are  derived  from  powers  of  W .  Power  queries  have 
been  successfully  used  for  a  number  of  problems,  see  again  [22],  including  the  phase  estima¬ 
tion  problem  that  will  be  discussed  in  the  next  section.  The  phase  estimation  algorithm, 
see  [8,  22],  is  at  the  core  of  many  quantum  algorithms.  It  plays  a  central  role  in  the  fast 
quantum  algorithms  for  factoring  and  discrete  logarithms  of  Shor  [27].  We  stress  that  for 
Shor’s  algorithm,  power  queries  can  be  implemented  by  a  number  of  elementary  quantum 
gates  that  is  polylog  in  n.  The  phase  estimation  algorithm  approximates  an  eigenvalue  of  a 
unitary  operator  W  using  the  corresponding  eigenvector,  or  its  approximation,  as  part  of  the 
initial  state.  The  powers  of  W  are  defined  by  Pi  =  2Z~1.  Therefore,  phase  estimation  uses 
queries  with  W\  =  W,  W2  =  W2,  W3  =  W2,2 ,  . . . ,  Wm  =  W2"'  \  It  is  typically  assumed, 
see  [8],  that  we  do  not  explicitly  know  W  but  we  are  given  quantum  devices  that  perform 
controlled- IF,  controlled- IF2,  controlled-IF2-',  and  so  on. 

For  the  Sturm-Liouville  eigenvalue  problem,  as  well  as  for  problems  studied  in  [25],  we 
will  use  the  matrix 


IF  =  exp  (iqMg)  with  i  =  \f—l  and  a  positive  7,  (31) 

where  the  n  x  n  matrix  Mq  was  introduced  in  Section  3.2  as  a  discretization  of  the  differential 
operator  Lg.  The  matrix  IF  is  unitary  since  Mq  is  symmetric. 

For  the  IFj  with  the  matrix  IF  of  (31)  we  modify  the  query  definition  in  equation  (21) 
and  assume,  as  in  [22,  Ch.  5],  that  for  each  j  the  Wj  is  one  quantum  query.  Accordingly, 
for  algorithms  that  can  be  expressed  in  the  form  (29),  the  number  of  power  queries  is  m, 
independently  of  the  powers  pj.  By  analogy  with  (22),  we  denote  their  error  by  equant(A,  m). 

Allowing  quantum  algorithms  of  the  form  (29)  with  power  queries,  we  define  the  power 
query  complexity  npower_query(e)  to  be  the  minimal  number  of  power  queries  required  to 
approximate  the  Sturm-Liouville  eigenvalue  problem  with  error  e,  i.e., 

npower— query  =  rnin{m  ;  3  \  such  that  equant(A,m)  <  e}. 

The  cost  of  one  power  query  is  denoted  by  cpower.  The  total  complexity,  comppower_query(e), 
is  the  defined  as  the  minimal  cost  of  a  hybrid  algorithm  in  the  same  way  as  for  bit  queries. 

We  will  use  the  phase  estimation  algorithm  as  a  basic  module  for  approximating  the 
smallest  eigenvalue  A (q).  As  shown  by  Abrams  and  Lloyd  [1],  the  phase  estimation  algorithms 
can  also  be  used  if  a  good  approximation  of  the  eigenvector  corresponding  to  the  smallest 
eigenvalue  is  known.  Such  an  approximation  is  obtained  by  the  algorithm  of  Jaksch  and 
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Papageorgiou  [18].  Combining  these  algorithms,  we  obtain  the  quantum  algorithm  that 
computes  the  smallest  eigenvalue  with  error  £  and  probability  |  using  ©(log  e”1)  power 
queries,  and  ©(log  e~l)  qubits. 

For  the  sake  of  completeness,  we  review  the  phase  estimation  problem  and  algorithm, 
the  results  of  Abrams  and  Lloyd  and  the  results  of  Jaksch  and  Papageorgiou  in  the  next 
subsections. 

5.3  Phase  Estimation 

Consider  W  defined  by  (31)  with  7  =  i.e., 


The  eigenvalues  of  W  are  elAAM</)/2  with  A j(Mq)  being  the  eigenvalues  of  the  n  x  n  matrix 
Mq  and  n  is  assumed  to  be  a  power  of  two.  These  eigenvalues  can  be  written  as  e2m<Pj ,  where 

<Pj  =  VAMq)  = 

are  called  phases.  We  are  interested  in  estimating  the  smallest  phase  tpi  (Mq),  which  belongs 
to  (0,1)  since  Ai (Mq)  €  [7r2, 7r2  +  1],  For  convenience,  we  renumber  and  normalize  the 
eigenvectors  of  Mq,  and  also  of  W,  as 

I  Vj)  =  Vn\zj+i(Mq )), 

for  j  —  0, 1, . . . ,  n  —  1.  We  will  use  { 1 2/y ) }  as  the  orthonormal  basis  of  the  space. 

Phase  estimation,  see  [22,  Section  5.2],  is  a  quantum  algorithm  that  approximates  the 
phase  (pi(Mq).  Note  that  to  compute  an  ^-approximation  of  Ai (Mq),  it  is  enough  to  com¬ 
pute  an  £ / (47r)-approximation  of  ipi(Mq).  The  original  phase  estimation  algorithm  has  been 
derived  for  the  initial  state  |O0m)l?/o),  where  m  is  related  to  the  accuracy  and  will  be  deter¬ 
mined  later,  and  | y0)  =  | y0(Mq))  is  the  eigenvector  of  the  matrix  Mq  corresponding  to  the 
smallest  eigenvalue  Ai  (Mq).  Abrams  and  Lloyd  [1]  showed  that  phase  estimation  can  still  be 
used  if  the  eigenvector  |r/0)  is  replaced  by  a  good  approximation  |^0)  as  the  initial  state. 

More  precisely,  expanding  |^0)  in  the  basis  of  the  eigenvectors  | yj),  the  initial  state  takes 
the  form 

n—  1 

|0>®”|<M  =  \ormY.dM- 

3= 0 

Using  m  Hadamard  gates,  we  place  the  first  register  in  an  equal  superposition,  which  gives 
the  state 

11  1  71—1 

1^)  =  ■"  \xi)\x^---\xrn)J2di\ yj)- 

'  11=0  12=0  Xm  =  0  3=0 

We  now  apply  the  controlled  quantum  gates,  see  (30),  to  create  the  state 

| ^2)  =  1 W2™'— 2  *  *  *  IVh) 

^  n—  1 

=  TTpE'tl’ft)  lw> 

V  3=0 
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with 


| rjj)  =  (jO)  +  e2^\l)^j  ®  (jO)  +  e27ri2^'|l)^  <g>  •  •  •  ®  (jO)  +  e2^i2m“Vi|1) 


i  1 


J]  J]  •  •  •  e27ri(a:i20+a:22l+-a:m2m'1)^|x1)|x2)  •  •  •  |xm) 


xi=0a;2=0  a;m=0 

2m_i 


1 


n—  1  /2m  — 1 


^=o 

see  also  [22,  p.  222],  Hence, 

i«  =  (  E  1 1»>- 

v  z  j= o  \  e=o  J 

The  inverse  Fourier  transform  performed  on  the  first  register  creates  the  state 

n— 1  /2m-l  \ 

^dj  (  S&j’t) \l)  1%)’ 

j= o  V  £=o  / 


where 


9(<PjJ)  = 


sin(7r(2myj-Q)e,rlwJ 


(¥>j~£2-m)(2m-l) 


1, 


2m  sin(7r((/?:7-— ^2— m)) 


if  ^  ^  2~”Y, 

if  ifj  =  2~m£. 


A  measurement  of  the  first  register  produces  the  outcome  j  with  probability 


n—  1 


pj  =  K \2\g{<pe,j)\2, 


£=0 


and  the  second  register  collapses  to  the  state 


n—  1 


E 


e=o 


VP] 


The  quantity 

A(0o,0i)  =  min{|x  +  fa  -  </>0\}  for  fafaeR 

#EZ 

is  dehned  in  [6]  and  is  the  fractional  part  of  the  distance  between  two  phases  (j) 0  and  <f>i .  It  is 
used  to  derive  the  relationship  between  the  approximation  error  and  the  success  probability. 
A  measurement  of  the  first  register  produces  an  outcome  from  the  set 

ft  =  {j:A(j/2l’,MM,))<k/2"'}, 

where  k  >  1,  with  probability 


n—  1 

P  r(<7fc)  =  ^2^2\deg((pe,j)  |2  >  |d|2  ^  |^(^i(M,),  j)|2  >  |d|2 

i&Qk  ^=o  jeQk 


\d\2 

2{k  —  1)  ’ 
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where  d  =  (yolV’o)-  For  k  —  1  the  probability  that 

A(j/2m,(pi(Mq))  <  2~m  is  bounded  from  below  by  -%\d\ 2 .  (32) 

7T 2 

The  proof  of  the  probability  bounds  can  be  found  in  [6,  22],  Using  this  fact,  the  authors 
of  [1]  conclude  that  as  long  as  |d|2  is  large  enough  or,  equivalently,  |^0)  is  close  enough  to 
\y0)  then  phase  estimation  can  be  used  to  approximate  the  phase  (pi(Mq)  with  probability 
close  to  8 / 7T2  =  0.81 .... 

We  stress  that  the  phase  estimation  algorithm  uses  m  power  queries.  In  addition  to  the 
cost  of  the  queries  there  is  a  quantum  operations  cost  proportional  to  at  most  m2,  which  is 
an  upper  bound  on  the  cost  of  the  quantum  inverse  Fourier  transform,  see  [22,  Section  5.2], 

5.4  Eigenvalue  and  Eigenvector  Approximation 

The  results  of  Jaksch  and  Papageorgiou  [18]  can  be  applied  to  efficiently  construct  a  good 
approximate  eigenvector  when  W  =  e 2 Mq  as  in  the  previous  subsection. 

The  matrix  Mq  =  Mqn^  has  been  derived  from  the  discretization  of  the  operator  hq  with 
mesh  size  hn  =  (n  +  l)-1.  Its  eigenvectors  are  also  eigenvectors  of  W  =  W^n\  and  we 
denote  them  here  by  | yf^),  where  j  =  0, 1, . . . ,  n  —  1.  We  want  to  approximate  Ai  (Mql'>)  = 
47T<^i(Mgn^)  but  we  do  not  know  the  corresponding  eigenvector 

l!/'”1)  :=  Ik"’)- 

The  expansion  of  \y^)  in  the  computational  basis  is  denoted  by 

n—  1 

l!/<n,>  =  Ei/Pli),  (33) 

3= 0 

/  ,  \  1/2 

Recall  that  uq  is  the  normalized,  ||uq||L2  —  (f0  u2(x)dx)  =  1,  eigenfunction  of  the  dif¬ 
ferential  operator  hq  that  corresponds  to  A (q),  and  uq  as  well  as  u'q  and  uq  are  uniformly 
bounded,  i.e.,  ||w9||oo,  1 1 1 1 oo  and  ||Ug||oo  are  0(1). 

Let  [f/O))  =  uq ( ( j  +  l)hn)\j)  be  the  vector  obtained  by  sampling  uq  at  the  dis¬ 

cretization  points.  Then  it  is  known,  see  [12,  19]  as  well  as  Remark  4.1,  that 

I'^-iiftlL  =  a,ld  (34) 

|a(9)«a1(m<“,)|  =  oul). 

Consider  a  coarse  discretization  of  hq  with  mesh  size  hno  =  (n o  +  l)^1  with  no  being  a  power 
of  two.  Assume  that 

no  — 1 

|||J(”»)>||2  =  1, 

3=0 
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approximates  the  eigenvector  \y^n°^)  that  corresponds  to  the  smallest  eigenvalue  of  the  matrix 
Mq"°')  such  that, 

II  |5(no)>  -  I y(no))  ||2  =  O(n0-2).  (35) 

We  place  the  vector  |ihn°'))  in  a  log  n0  qubit  register.  As  explained  in  Section  4.3,  we  can 
compute  |5(no'))  on  a  classical  computer  with  cost  of  order  Uq. 

For  n  =  2sno,  we  construct  an  approximation  | z^)  of  | y^)  by  hrst  appending  s  qnbits, 
all  in  the  state  |0),  to  |5^n°^)  and  then  performing  the  Hadamard  transformation  on  each  one 
of  these  s  qnbits,  i.e., 


where  ’s  denote  the  coordinates  of  |5^n°^)  in  the  computational  basis,  and  g(j )  =  |_j /2SJ . 
The  effect  of  g  is  to  replicate  2s  times  the  coordinates  of  |5^n°^).  As  in  Jaksch  and  Papageor- 
giou  [18],  we  use  the  vector  | z^)  as  part  of  the  input  to  the  phase  estimation  algorithm. 

Let  d =  {y^\z^}.  We  show  that  \d^\2  can  be  made  arbitrarily  close  to  one  by 
choosing  a  sufficiently  large  Uq.  Hence,  we  can  make  the  success  probability  of  the  phase 
estimation  algorithm  at  least  equal  to  |. 

Consider  two  different  expansions  of  | z^), 

n—  1 

y”,>  =  EaPU>  (37) 

3=0 
n—  1 

I  *<”>)  =  EAkb).  (38) 

i-'i 

The  hrst  expansion  is  in  the  computational  basis  (|j)}  and,  by  (36), 

u for  3  =  0, 1,  •  •  • ,  n  -  1, 

while  the  second  expansion  is  with  respect  to  the  eigenvectors  of  Mq'} .  Note  that  =  d^ 
and  clearly  =  1-  Equation  (38)  implies 


z(n))  ~  | y{n))  =  (d(n) 


n—  1 

Div^+EWt 


3= 1 


(39) 


Taking  norms  on  both  sides  we  obtain 

n—  1  n—  1 

II  |j/”))  -  |z(”>)  \  \l  =  |rf<"»-l|2  +  ElA,|2  >  EP”1!2  =  1-M<n)l2-  («) 

3= 1  3= 1 
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We  now  bound  the  left  hand  side  of  (40)  from  above.  Using  the  expression  (33)  for  | 
and  the  definition  of  | z^),  see  (36),  (37),  we  have 


n—  1 


|||!/(”,)-|«(”»)||2  = 

3=0 
n—  1 

-  E 


3=0 


uq(U  jl  4~  l)hno)  .  (n)  _  ^g(j) 


III  u^)b 


V¥\\\u(no) 


where,  by  (34)  and  (35),  we  have 

n—  1 

E  iAfi2  =  °(;4) 

j=0 

Applying  the  triangle  inequality,  we  get 

/  n—  1 

||i»(”)>-i^<”)>||2  <  E 


and 


n—  1 

2*0(0 

3=0 


<3=0 


uq((j  +  l)hn)  uq{{g(j)  +  l)h 


noV 


ll|W)lh 


v^*|||£7(n°) 


1/2 


+  o(h: 


no> 


(41) 


The  definition  of  |f/(-n'))  and  the  fact  that  the  derivative  of  uq  is  Lipschitz3  with  the  uniform 
Lipschitz  constant  imply  that  |||C/^)||2  =  \/n(l  +  0(hn)).  Hence,  the  square  of  the  term  in 
the  parentheses  above  is  equal  to 

n—  1 


l 


n 


I Uq((j  +  l)^n)(l  +  0(hn ))  —  Uq{{g(j)  +  l)/lno)(l  +  O(hno))\2. 


(42) 


3=0 


Since  uq  is  continuous  with  a  bounded  first  derivative,  we  have  that 

Uq{x2,j)  =  Uq(x ij)  +  0( \x2J  -  Xl,j\),  (43) 

where  x2,j  =  ( j  +  1  )hn  and  xhj  =  (g(j)  +  1  )hno,  j  =  0, 1, . . . ,  n  -  1.  Let  j_j'/2sJ  =  j/2s  -  a 
with  a  G  (0, 1).  Then 

j  +  1  J'/2S  +  1  -  a 


1*2*2,/  *2*1,/ 1 


2sn0  T  1  no  +  1 

=  J  (2°n0  +  l)2s(n0  +  1)  +  °(/lno)  =  °^n°' 
Using  (42),  (43)  and  the  triangle  inequality,  we  obtain  from  (41)  that 

ili</’*,>-h“)}||2  =  o(hm)  < 


n0  +  1 

for  some  positive  number  c  independent  of  n  and  no-  Combining  this  with  (40)  we  finally 
conclude  that 


|d(n)|2  >  1  - 


(no  +  l)2 

and  dSn^  can  be  made  arbitrarily  close  to  one  by  taking  a  sufficiently  large  no- 


(44) 


3 A  function  /  :  [0,1] 
all  x,  y  €  [0, 1], 


is  Lipschitz  if  there  is  a  number  L  >  0  such  that  | f(x)  —  f(y)\  <  L\x  —  y\  for 


5.5  Quantum  Algorithm  for  the  Smallest  Eigenvalue 

We  combine  the  results  of  the  previous  two  subsections  to  derive  a  quantum  algorithm  for 
computing  an  e- approximation  of  the  smallest  eigenvalue  with  probability 

We  choose  the  parameters  for  the  phase  estimation  algorithm.  Without  loss  of  generality 
we  assume  that  e_1  is  an  even  power  of  2,  that  is  e"1  =  2m  with  an  even  m.  We  set 
n  =  £-1/2  =  2m/2  and  we  will  be  working  with  the  matrix  Mq'] .  The  index  no  =  2k°  is 
chosen  as  the  smallest  power  of  two  for  which 


8  /  c2  \  3 

it2  \  (n0  +  l)2/  _  4’ 


(45) 


where  the  number  c  is  from  (44).  Clearly,  Uq  =  0(1).  Without  loss  of  generality  we  assume 
that  \m  >  k0  =  log  n0,  i.e.,  we  assume  that  £  is  sufficiently  small.  We  finally  set  s  =  | m  —  k0 . 
We  then  compute  2-n°))  on  a  classical  computer  as  in  Section  4.3  with  cost  0(1)  function 
values  and  operations. 

i  »  ,(n) 

We  run  the  phase  estimation  algorithm  for  the  matrix  W  =  e^Mq  with  the  initial  state, 
see  (36), 

|0)(g>m|  ~(n))  =  |0)®m|~(no)^  . 

Let  j  be  the  outcome  of  the  phase  estimation  algorithm.  We  finally  compute 


Xj  =  4vr  j  2~m 

as  an  approximation  of  the  smallest  eigenvalue  A (q).  We  have 

A,  -  \(q)  =  y  -  A ,(MW)  +  A ,(MW)  -  A (q) 


=  to  1  -L  -  +  0(e). 


From  (32)  we  know  that 


i  -  <Pi (MW)  )  <  £  with  probability  —  \d^\2. 


By  (44)  and  the  definition  of  no  we  have 

4l d{n)\2  > 

Hence, 


7T- 


n2  V  (n0  +  l)2 


TTC 


>4 

-  4 


\X3  —  A(g)|  =  0(e)  with  probability  at  least 


The  computation  of  A  j  requires 


m  +  k0  +  s  =  \  rn  =  f  log  e 


-l 
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qubits,  m,  =  log  e_1  power  queries,  plus  a  number  of  quantum  operations  proportional  to 
m 2  =  log2  1 .  This  yields  nP°wer-query(+)  —  0(log  £~l).  A  lower  bound  on  npower"query(e)  of 
the  same  order  is  proved  in  [5].  Hence, 

^power-query (^ )  =  Q(1()g  £~ly 

We  summarize  the  results  of  this  section  in  the  following  theorem. 

Theorem  5.3.  The  Sturm- Liouville  eigenvalue  problem  can  be  solved  in  the  quantum  setting 
with  power  queries  by  the  phase  estimation  algorithm  applied  to  the  discretized  matrix  of  the 
differential  operator  Lq  with  the  initial  state  given  as  an  approximate  eigenvector  computed  by 
the  Jaksch  and  Papageorgiou  algorithm.  This  quantum  algorithm  approximates  the  smallest 
eigenvalue  A (q)  with  error  e  and  probability  |  using 

•  |  log  e-1  +  0(1)  power  queries, 

•  0(1)  function  values  and  classical  operations, 

•  0(log2e_1)  quantum  operations  besides  the  power  queries,  and 

•  |  logs-1  +0(1)  qubits. 

Furthermore, 

^power-query  =  QQog  £-l^ 

and 

^(Cpower  log  £^)  =  C0mppOWer‘qUery (e)  =  O  (cpower  log  £_1  +  C  +  log2^1)  . 

5.6  Qubit  Complexity 

In  this  section  we  address  the  qubit  complexity,  compqub(e),  which  is  defined  as  the  minimal 
number  of  qubits  required  to  approximate  the  smallest  eigenvalue  with  error  £  and  prob¬ 
ability  |  by  quantum  algorithms  of  the  form  (21).  Clearly,  compqub(e)  is  upper  bounded 
by  flog  e^1  +  0(1)  since  that  many  qubits  are  used  by  the  phase  estimation  algorithm  of 
Section  5.5.  Observe  that  the  cost  of  the  classical  algorithm  computing  \z^no'>)  as  well  as  its 
quantum  simulation  [22,  p.  189-193]  is  constant  since  n0  is  bounded  by  a  constant  due  to 
(45). 

We  turn  to  a  lower  bound  on  compqub(e).  Based  on  the  results  obtained  in  this  paper,  it 
is  easy  to  see  that  the  number  of  qubits  necessary  to  solve  our  problem  must  be  proportional 
at  least  to  roughly  lloge-1.  Indeed,  assume  that  there  is  a  quantum  algorithm  of  the 
form  (21)  that  computes  A (q)  with  error  £  and  probability  |,  and  uses  k(e)  qubits.  This 
algorithm  can  use  arbitrary  quantum  queries,  assuming  that  each  quantum  query  is  based  on 
at  most  2k ^  function  evaluations  of  q.  Note  that  this  holds  for  bit  queries,  as  well  as  for  the 
power  queries  studied  in  this  paper.  Then  such  an  algorithm  can  be  simulated  by  a  classical 
algorithm  that  uses  at  most  2k &  function  evaluations  of  q.  From  Theorem  3.2  we  know  that 
2fcO  =  f2(e-1/2)  and  therefore  k(e)  >  1  log  e~l  +  f2(l).  Hence,  the  qubit  complexity  is  lower 
bounded  by  \  log  e~l  +  0(1) .  This  proves  the  following  theorem. 
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Theorem  5.4.  The  qubit  complexity  of  the  Sturm- Liouville  eigenvalue  problem  in  the  quan¬ 
tum  setting  with  bit  or  power  queries  is  bounded  by 


|  log  £  1  +  0(1)  <  compqub(e)  <  f  log  £  1  +  0(1). 
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